/*matvec.c
   Vector and matrix allocation functions.
   Matrix initialization function.
*/
#include <stdio.h>
#include <stdlib.h>
#include "matvec.h"

/* allocate memory for a vector of int of length n*/
int *ivector(size_t n)
{
   int *a;
   a = (int *)malloc(n * sizeof(int));
   return a;
}

/* allocate memory for a vector of double of length n*/
double *dvector(size_t n)
{
   double *a;
   a = (double *)malloc(n * sizeof(double));
   return a;
}

/* realloc memory for a vectoy of double of length n*/
double *dvector_realloc(double *a, size_t n)
{
   a = (double *)realloc(a, n * sizeof(double));
   return a;
}

/* allocate memory for a vector of float of length n*/
float *fvector(size_t n)
{
   float *a;
   a = malloc(n * sizeof *a);
   return a;
}

/* allocate memory for a matrix of int of size rows by columns*/
int **imatrix(size_t rows, size_t columns)
{
  int **a;
  size_t i,j;
  a = malloc((rows+1) * sizeof *a);
  if (a==NULL)
      return a;

  for(i=0; i<rows; i++) {
    a[i] = ivector(columns);
    if (a[i] == NULL) {
       for(j=0; j<i; j++)
         free_ivector(a[j]);
         free(a);
         return NULL;
    }
  }

  a[rows] = NULL;
  return a;
}

/* allocate memory for a matrix of double of size rows by columns */
double **mydmatrix(size_t rows, size_t columns)
{
  double **a;
  size_t i;
  size_t j;

  a = malloc((rows+1) * sizeof *a);

  if (a==NULL)
      return a;

  for(i=0; i<rows; i++) {
    a[i] = dvector(columns);
    if (a[i] == NULL){
       for(j=0; j<i; j++)
        free_dvector(a[j]);
        free(a);
        return NULL; 
    }
  }

  a[rows] = NULL;

  return a;
}

/* reallocate memory for a matrix of double of size rows by columns */
double **mydmatrix_realloc(double **a, size_t rows, size_t columns)
{
  size_t i;
  size_t j;

  a = (double**)realloc(a, (rows) * sizeof(double*));

  if (a==NULL)
      return a;

  for(i=0;i<rows; i++)
	a[i] = NULL;

  for(i=0; i<rows; i++) {
//    a[i] = dvector_realloc(a[i], columns);
	a[i] = (double *)realloc(a[i], columns*sizeof(double));
    if (a[i] == NULL){
       for(j=0; j<i; j++)
        free_dvector(a[j]);
        free(a);
        return NULL;
    }

  }

//  a[rows] = NULL;

  return a;
}

/* allocate memory for a matrix of float of size rows by columns*/
float **fmatrix(size_t rows, size_t columns)
{
  float **a;
  size_t i;
  size_t j;

  a = malloc((rows+1) * sizeof *a);

  if (a==NULL)
      return a;

  for(i=0; i<rows; i++) {
    a[i] = fvector(columns);
    if (a[i] == NULL){
       for(j=0; j<i; j++)
        free_fvector(a[j]);
        free(a);
        return NULL; 
    }
  }

  a[rows] = NULL;

  return a;
}

/*free memory allocated by ivector()*/
void free_ivector(int *a)
{
	free(a);
}

/*free memory allocated by dvector()*/
void free_dvector(double *a)
{
	free(a);
}

/*free memory allocated by fvector()*/
void free_fvector(float *a)
{
	free(a);
}

/*free memory allocated by imatrix()*/
void free_imatrix(int **a)
{
	size_t i;
	 
	if (a==NULL)
		return;
	
	for (i=0; a[i]!=NULL; i++)
		free (a[i]);

	free(a);
}

/*free memory allocated by mydmatrix()*/
void free_mydmatrix(double **a)
{
	size_t i;
	
	if (a==NULL)
		return;
	
	for (i=0; a[i]!=NULL; i++)
		free (a[i]);

	free(a);
}

/*free memory allocated by mydmatrix()*/
void free_fmatrix(float **a)
{
	size_t i;
 	
	if(a==NULL)
		return;
	
	for (i=0; a[i]!=NULL; i++)
		free (a[i]);

	free(a);
}

/* initialize mtarix of type double with 0 */
double **initial_matrix(int m, int n)
{	
	double **a;
	size_t i, j;

	a = mydmatrix(m,n);
	
	if(a == NULL){
		fprintf(stderr, "out of memory in function "
                             "%s, file %s, line %d\n",
                       "initial_matrix()", __FILE__, __LINE__);
		return NULL;
	}
	
	for (i=0; i<m; i++)
		for(j=0; j<n; j++)
			a[i][j] = 0.0/*i*n+j*/;
	
	return a;
}
	
double *initial_vector(int n)
{
	double *a; 
	int i;

	a = dvector(n);
	if(a == NULL){
		fprintf(stderr, "out of memory in function "
                             "%s, file %s, line %d\n",
                       "initial_vector()", __FILE__, __LINE__);
		return NULL;
	}
	for(i=0; i<n; i++)
		a[i] = 0.0;

	return a;

}


